function sigma = Zhang(img, k)
    [I_rows, I_cols] = size(img);
    img = double(img);
    % 计算图像原本的p值
    img(isnan(img)) = 0;
    ftf = img'*img;
    [~,S,~] = svd(ftf);
    p = trace(sqrt(S))/I_cols;
    % 计算图像添加辅助噪声的p值
    
    sigma_add = 20;
    img_add = img + sigma_add*randn(I_rows, I_cols);
    ftf = img_add'*img_add;
    [~,S,~] = svd(ftf);
    p_add = trace(sqrt(S))/I_cols;
    H = (p_add - p)/k;
    sigma = (sigma_add^2 - H^2)/(2*H);
end
